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ABSTRACT 

We construct models of the structural evolution of super-Earth- and mini-Neptune-type exoplanets 
with H 2 -He envelopes, incorporating radiative cooling and XUV-driven mass loss. We conduct a pa¬ 
rameter study of these models, focusing on initial mass, radius, and envelope mass fractions, as well 
as orbital distance, metallicity, and the specific prescription for mass loss. From these calculations, 
we investigate how the observed masses and radii of exoplanets today relate to the distribution of 
their initial conditions. Orbital distance and initial envelope mass fraction are the most important 
factors determining planetary evolution, particular radius evolution. Initial mass also becomes im¬ 
portant below a “turnoff mass,” which varies with orbital distance, with mass-radius curves being 
approximately flat for higher masses. Initial radius is the least important parameter we study, with 
very little difference between the hot start and cold start limits after an age of 100 Myr. Model sets 
with no mass loss fail to produce results consistent with observations, but a plausible range of mass 
loss scenarios is allowed. In addition, we present scenarios for the formation of the Kepler-11 planets. 

Our best fit to observations Kepler-llb and Kepler-llc involves formation beyond the snow line, after 
which they moved inward, circularized, and underwent a reduced degree mass loss. 


1. INTRODUCTION 

Recent exoplanet surveys, particularly the Kepler Mis¬ 
sion, have discovered hundreds of confirmed exoplan¬ 
ets and thousands of candidates, prompting studies of 
the general structures and evolution of planetary bodies. 
Notably, the most common class of planets discovered 
by these surveys are “super-Earths” or “mini-Neptunes” 
with radii of ~l-4 R®, with an apparent modal value of 
~2-3 Rrn, for w hi ch the re is no analog in our Solar Sys¬ 
tem dBurke et all l2014h . This appears to be true even 
for a sample set of planets corrected for completeness 
(|Silburt et al.lf20T4l ). For objects for which masses have 
been measured, the inferred densities of planets with 
radii of 2-4 R® are too low to be accounted for without 
the presence of a voluminous hydrogen-helium envelope. 
However, their densities fall in a wide range so that plan¬ 
ets of the same mass may have very different radii; they 
may be terrestrial, or they may have envelopes compris¬ 
ing up to tens of p ercent of their mass (|Marcv et al.ll20l4 
iHowe et al.ll2014r) . 

A small number of very-low-density planets have 
been observed with radii >5 R®, but masses 

less than that of Neptune, such as K epler-18c 
and -18d (|Cochran et al.| 1 2 01 ID. Kepler-51b (IMasudal 
120141) Kepler-79d (IJontof-Hutter et all I2014D. Kepler- 
87b (|Ofir et all 120141) . and Kepler-89e (|Masuda et al.1 
I2013D . as demonstrated variously by radial velocity 
(RV) and transit timing variation (TTV) measurements. 
These objects have so far been found only at larger or¬ 
bital distances of 0.1-0.2 AU, which, along with the¬ 


oretical considerations, suggests that they typically do 
not retain a H 2 -He atmosphere at <0.1 AU due to rapid 
mass loss or other factors. 

Many of these planets lie within 0.2 AU of their host 
stars, suggesting at the very least that they may have 
been significantly sculpted by stellar flux-driven mass 
loss. This is true even at larger masses up to the “hot 
Jupiters”; for example, HD 209458b appears to be ex¬ 
periencing measurable mass los s today on the order o f 
2 x 10 10 g s- 1 (10 -10 M® yr -1 ) (lEhrenreich et al.ll2008D . 
and some apparently-rocky planets may, in fact, have 
formed with significant envelopes that have since been 
entirely lost. Therefore, a self-consistent model of extra¬ 
solar planets with hydrogen-helium atmospheres should 
incorporate the effects of Kelvin-Helmholtz cooling and 
mass loss to determine whether a particular planet has 
likely been subjected to significant alteration by mass 
loss. 

Structural models of planets indicate that a planet’s 
radius may serve as a reasonably good proxy for esti- 
mati ng its composition—eith er the envelope mass frac¬ 
tion dLopez fe Fortnevl I2014D or the mass of the enve¬ 
lope itself (|Howe et afll2014D . However, there are many 
degeneracies involved in these estimates, and uncertain¬ 
ties in the data cause further difficulties in interpreta¬ 
tion. However, combining good mass, radius, and age 
measurements with evolutionary models may allow less 
ambiguous interpretations to be made. 

An evolutionary history of a planet incorporating the 
effects of radiative cooling and mass loss (particularly 
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mass loss due to XUV ph otoevaporatio n) could also give 
clues to planet formation. iRibas et~ahl (12005T) found that 
solar-type stars produce the greatest amount of ultra¬ 
violet radiation when they are young, declining rapidly 
after 100 Myr. Therefore, if such an evolutionary model 
can determine whether a given planet could have sur¬ 
vived in its present location during the high-mass-loss 
early phases of its existence, it would provide further in¬ 
formation about where and when planets form in their 
proto-planetary disks. 

There is some uncertainty in the effect of mass loss. 
iLopez et al] ( 201 2;) constructed evolutionary models in¬ 
corporating cooling and mass loss for the specific cases 
of the planets of the Kepler-11 system. They specifi¬ 
cally adopted an energy-limited XUV-driven mass loss 
prescription, finding that the five nearest planets to the 
star (within 0.25 AU) have undergone significant mass 
loss, but that only the inne rmost has potentially los t its 
entire envelope. However, ILopez &; Fortnevl (12013D . in 
a similar study of the Kepler-36 system, found that the 
transition in orbital distance from retaining nearly all of 
the envelope to losing all or nearly all of it is relatively 
sharp for a give n core mass ._ 

Previously, in lHowe et al.l ( 2014 1. we studied the struc¬ 
tures of exoplanets with gaseous envelopes, finding that 
observed extrasolar planets may have a wide range of 
core masses and envelope masses with little correlation 
between the two. Therefore, a wide range of models is 
needed to fully characterize the exoplanet population. 
For evolutionary models, this problem is compounded 
because the planets may have a wide range of orbital 
distances and orbit stars with a range of spectral types, 
resulting in a wide range of temperatures of their atmo¬ 
spheres. Also, the initial radii of the planets are not cer¬ 
tain, although_somewoiThns been done to mo del them, 
e.g. IMordasini et al.l (|2012aD : I.Tin et al.l (|2014h . Finally, 
the detailed properties of the atmospheres, such as cloud 
cover, redistribution of heat, and the processes that con¬ 
tribute to mass loss, are poorly understood. 

While there is a rich history of exoplanet structure 
modeling and a growing body of work on planetary evo¬ 
lution, much remains to be done. Uncertainties must be 
addressed and plausible parameters set to produce use¬ 
ful models. Ideally, a range of models would be produced 
with different parameters, but in the case of atmospheric 
properties, which are not fully understood, a simple pa¬ 
rameterization is perforce necessary. Also, these tech¬ 
niques to study evolution have only r ecently been applie d 
to super-Earths and mini-Neptunes (|Lopez et al.ll201§ ). 

Our paper investigates the evolution of theoretical 
super-Earth and mini-Neptune models in different plan¬ 
etary systems over billions of years. We focus particu¬ 
larly on finding the boundaries of the parameter space 
in which a planet’s envelope is lost completely and on 
estimating the initial conditions of observed planets. 

Section [2] summarizes the previous work on this topic. 
Section [3] describes our planetary model. We demon¬ 
strate our code in Section [4] by providing internal den¬ 
sity and pressure profiles over time for example models. 
Then, in Section 0 we show the resultant mass loss rates 
and investigate the effects of varying the prescription for 
mass loss. In Section O we conduct a parameter study of 
the most important physical parameters that go into the 
calculations, constructing models with a range of initial 


masses, radii, envelope mass fractions, orbital distances, 
metallicities, and mass loss prescriptions. In Section [71 
we attempt to fit the evolutionary history of the plan¬ 
ets in the Kepler-11 system. Section [8] provides example 
analyses of the evolution of transit spectra of our models, 
and we present general conclusions in Section [9] 

2. PREVIOUS WORK ON EVOLUTIONARY MODELS AND 
ATMOSPHERIC ESCAPE 

One importa nt early study o f atmo spheric escape was 
undertaken by iHubbard et al.1 (|2007f) . T hey estimated 
mass l oss by scaling the XUV flux model of Watson et al. 
(1981) to apply to highly-irradiated giant planets. They 
applied the tidal effect to b oth the Watson models and 
the model of iBaraffe et al.l (12004 1. They also added an 
efficiency parameter—the fraction of the total stellar flux 
that contributes to energy-limited escape, which is less 
than the total XUV flux due to losses through dissocia¬ 
tion and ionization, in contrast with the Baraffe model, 
which used the full XU V flux. F urt her, they applied 
the XUV flux model of IRibas et all ( 2005 ). The over¬ 
all efficiency parameter (the fraction of total stellar flux 
that goes into the evaporation of the atmosphere) at 4.5 
Gyr is ~ 10 -4 for the Baraffe model, compared with 
~ 1CU 6 for the Watson model. The efficiency parameter 
is proportional to the XUV flux, which is considered to 
decrease at a rate proport i onal t o t -123 after 100 Myr. 
Therefore, IHubbard et all (120071) found a mass l oss rate 
a factor of 100 smaller than IBaraffe et ahl (|2004f l . which 
eliminates the problem of improbable model parameters 
for HD 209458b found by that study. They found that 
mass loss becomes significant for highly-irradiated plan¬ 
ets in general only for planets of about the mass of Saturn 
and sma ller. _ 

IMordasini et al.1 (|2012aD and IMordasini et al.l (|2012bf ) 
modeled a wide range of planets with radioactive heating 
and Kelvin-Helmholtz cooling, but not mass loss. No¬ 
tably, they also modeled planet formation with a core¬ 
accretion model incorporating migration, providing esti¬ 
mates for initial radii at the end of the accretion phase, 
and modeled radiative cooling with gray atmospheric 
boundary conditions in the Eddington approximation. 
However, they set a minimum orbital distance of 0.1 
AU (as opposed to an observed limit of about 0.015 AU 
for 55 Cancri e) so that the strongest irradiation effects 
are avoided. Using this framework, they model the for¬ 
mation and evolution of planets with masses from 1 to 
10,000 M® (0.003 to 30 Mj ) and envelope mass frac¬ 
tions ( fenv) from < 1% to > 99%. These objects have 
radii of < 20R® (< 2 Rj) once the accretion phase ends 
and follow a relatively narrow mass-radius relation after 
100 Myr for which radius is a good first-order proxy for 
the fraction of mass in the envelope, f env and mass and 
radius together are sufficient to compute f env with rela¬ 
tively high precision. However, in all cases, they assumed 
a low-entropy “cold start” scenario in which the energy 
of accretion is thermally radiated on a timescale that is 
short compared with the accretion timescale, so that the 
effective temperature is always relatively low. This is in 
contrast to a “hot start,” in which the energy of accretion 
is not radiated efficiently, resulting in a much greater en¬ 
velope entropy and effective temperature when accretion 
ends. 

Meanwhile, iMurrav-Clav et ahl (120091) modeled mass 
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loss rates for hot Jupiters as a function of XUV flux, 
taking into account additional processes such as ioniza¬ 
tion, radiative losses, plasma recombination, tidal grav¬ 
ity, pressure confinement by the stellar wind, and mag¬ 
netic pressure. They found that mass loss rates for 
hot Jupiters are approximately energy-limited for lower 
fluxes (< 10 4 erg cm -2 s _1 or about the flux received by 
a planet 0.2 AU away from a young G dwarf at <100 Myr 
in the Watson model). The nearly energy-limited mass 
loss follows a relation of to oc F ( ^f Jv . However, at higher 
fluxes, this relation becomes flatter, roughly to oc F^fry 
d ue to rad iative and recombination losses. 

lOwen fc Wul (12013 1 examined a detailed model of mass 
loss going beyond the simple, energy-limited prescrip¬ 
tion. They incorporated X-ray emissions in addition 
to XUV, and they accounted for the recombination and 
line cooling processes that can reduce the mass loss effi¬ 
ciency. They modeled both Jupiter-mass and Neptune- 
mass planets based on a structure model using the MESA 
code with a range of initial radii from about 10 to more 
than 20 R®—notably including Neptune-mass models 
with large cores of 12-15 M®. With this mass loss pre¬ 
scription, they found significant mass loss approaching 
most of the hydrogen fraction (and a consummate de¬ 
crease in radius) for highly irradiated planets closer than 
-0.1 AU. 

iJin et all (120141) added a tmospheric escap e to the pop- 
ulation synthesi s model of iMordasini et all (|2012afl and 
IMordasini et all (l2012bD , using multiple mass loss pre¬ 
scriptions, including XUV-driven mass loss only, EUV 
plus X-ray-driven mass loss, purely energy-limited mass 
loss, and mass loss that is moderated by cooling radiation 
in the high-EUV regime, in which case the mass outflow 
is at the sound speed. In order to better study evapora¬ 
tion, they shifted the entire planet population inward by 
0. 04 AU, placi n g the inner edge at 0.06 AU. 

iLopez et all (|2012f ) modeled thermal evolution and 
XUV-driven mass loss for low-mass, super-Earth/mini- 
Neptune planets, rather than giant planets. In these 
models, a planetary structure is modeled with a fully 
convective envelope, with a structure determined by the 
specific entropy, which cools radiatively over time. They 
incorporated radioactive heating and the thermal iner¬ 
tia of a rocky core. The radiative flux is determined 
fr om a grid of model atmos pheres based on the models 
of iNettelmann et alJ (120111 4 . They used a “hot start” 
with a large initial entropy. The XUV-driven mass loss 
in this model is computed by an energy-limited formula 
with an efficiency of 10j)g°% applied to the XUV flux 
specifically and an XUV-photosphere at a pressure level 
from 0.1 nbar to 1 0 nbar. A more precise calculation 
from the iRibas et all (120051) re sults in an efficiency pa¬ 
rameter under the iHubbard et all (|2007f ) formulation of 
3. 3 x 10~ 6 at 4 5 Gy r- 

ILopez et all (|2012i l apply their models to the super- 
Earths and mini-Neptunes, Kepler-llb through f, in an 
attempt to determine their initial compositions. In their 
best fits, all five planets formed with hydrogen-helium en¬ 
velopes comprising at least a few percent of their masses, 
and they all subsequently lost at least 20% of the enve¬ 
lope masses after 100 Myr, with Kepler-llb being almost 
entirel y stripped. 

In ILopez fc Fortnevl (120131 1 , these same models were 


applied to Kepler-36b and -36c. Kepler-36b (mass 4.5 
M®) appears to be entirely rocky, while Kepler-36c (mass 
8.1 M®) has a density less than that of water, suggesting 
a voluminous hydrogen-helium atmosphere, despite be¬ 
ing in orbits at similar distances (an apparent 6:7 mean 
motion resonance). They model mass loss over 5 Gyr 
and find that the density dichotomy can be explained by 
a 65% more massive core for Kepler-36c, which results in 
a mass loss timescale —3 times as long as for Kepler-36b, 
al lowing it to retain an en velope over its lifetime. 

ILopez fc Fortnevl (|2014T) apply these models to the gen¬ 
eral case, producing models with a wide range of masses 
(1-20 M®), envelope mass fractions, f env (0.01%-20%), 
and levels of stellar irradiation (0.1-1000 times Earth’s ir¬ 
radiation, F®). They find that, to first order, planetary 
radius is dependent on only f env for larger envelopes of 
f env > 1%. They also propose a dividing line of 1.75 R® 
between planets that are entirely rocky and those with 
g aseous envelopes. 

I Wolfgang fc Lopezl (12014 1 apply the models of 
ILopez et al.l (1201214 and subsequent papers to the ob¬ 
served distribution of planetary radii from the Kepler 
mission (selected for completion) using a hierarchical 
Bayesian model. They find that the observed population 
has a distribution of f env with a mean of 0.7% and a 
standard deviation of 0.6 dex, corresponding to a modal 
r adius of —2 .0-2.5 R m . 

IWeiss fc Marcv (12014 ) fit a mass-radius function to 65 
exoplanets with radii less than 4 R®. In doing so, they set 
a boundary of 1.5 R® between rocky planets and those 
with gaseous envelopes. For rocky planets their best fit 
is a density of p = 2.43 + 3.39(i?/i?®) g cm -3 , while for 
planets with gaseous envelopes, their fit for the average 
radius is a mass of (M/M®) = 2.69(I?/R®) 0 - 93 . Their 
m ass, radius, and density measurements are taken from 
iMarcv et al.l (|2014f l. 

In our previous work. lHowe et all (|2014f) . we presented 
structural models of super-Earths and mini-Neptunes 
without thermal evolution, but at a range of envelope en¬ 
tropies, potentially representing different ages or levels of 
irradiation. In Figure Q] we plot mass-radius curves for 
a few selected compositions with observed mass and ra¬ 
dius measurements of planets. This figure includes mod¬ 
els with pure iron and terrestrial compositions plus ter¬ 
restrial compositions with H 2 -He envelopes comprising 
0.1%, 1%, and 10% of their masses with entropies of 6.0 
ks per baryon. We also found that radius is a proxy for 
composition, but proposed that it is more dependent on 
the envelope mass, m env . On the other hand, planets 
with similar masses may vary greatly in radius. Assum¬ 
ing an Earth-like core, the range of observed radii may 
indicate envelope masses from 0.01 M® to several M® 
for a given total mass between 2 and 10 M® (f env on the 
order of 0.1% to 50%). 

These results are updated and summarized in Figure O 
for the case of an Earth-like composition core and an 
envelope including a radiative layer. Figure [2] shows 
the relationship between envelope fraction and plane¬ 
tary radius for planets with observed masses and radii 
for envelope entropies of 5.5, 6.0, and 6.5 kb per baryon. 
These entropies a re re presentative of those found by 
iLonez fc Fortnevl (120144 . We consider “rocky” or ter- 

1 A core-mantle structure of 32.5% Fe and 67.5% MgSiC> 3 . 
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Fig. 1.— Mass-radius diagram for low-mass exoplanets with ob¬ 
served masses and radii, compared with mass-radius curves for se¬ 
lected terrestrial and gaseous composit ions with an envelope en¬ 
tropy of 6.0 kg/B. Based on data from IMowe et all <120 1 41 1 . 



Fig. 2.— Our calc ul ated envelope fractions (based on numerical 
data from IHowe et al~1l2014 *) for planets with observed masses and 
radii versus radius. We compute models with envelope entropies 
of 5.5, 6.0, and 6.5 fcj, per baryon. 


restrial planets to be those with a convective envelope 
with a mass fraction f env < 0 . 0001 , excluding the (much 
lower mass) radiative atmosphere. Often, a bare core 
will fit the observed radius. We find that the values are 
closely correlated, with radius increasing rapidly for high 
envelo pe fractions. This result is qualitatively similar to 
that of lWolfgang fe Loped (|2014f l , but with greater scat¬ 
ter and a higher apparent median value of f env of several 


percent as opposed to 0.7%. However, this observational 
sample is not necessarily complete. We also find that 
the lowest-density planets with radii larger than Nep¬ 
tune have f env > 10%. (Com pare Uranus and Neptune, 
w hich have = 1 0%-20%, (Helled et al.ll201lll 
Howe ct al. (2014) also presented fits with models with 
pure water ice cores, which produced similar results, ex¬ 
cept that the envelope masses were smaller by tens of 
percent; for small planets (R < 2 i?®), no envelope was 
needed to fit the data. Some additional uncertainty re¬ 
mained because a range of envelope entropies is possible. 
We found that a range of 5.5-6.5 fcs per baryon in en¬ 
tropy was associated with an uncertainty of tens of per¬ 
cent in envelope masses, roughly on the same order as 
the observational uncertainties. 


3. EVAPORATION AND EVOLUTION MODEL 


We compute pl anet a ry st ru cture profiles using the code 
we describe in IHowe et~a!1 (12014H . In order to model 
thermal evolution, we add to the equations of hydrostatic 
equilibrium, 


dP_ 

dr 

dm 

dr 


Gmp 


= 4-7T r 2 p, 


( 1 ) 

( 2 ) 


equations of thermal transport and energy production: 


dT 

dr 



T dP 
P dr 


dL 

dr 


ditr 2 pe, 


(3) 

(4) 


where T is the temperature, L is the luminosity, 7 is 
the adiabatic index, e is the energy production per unit 
mass, e.g. from radioactive decay. A convective thermal 
gradient in the atmosphere is assumed. 

We model radiative cooling by relating the energy bud¬ 
get to the change in envelope entropy: 


’ M C ore 


dt 


J radio 


C v Mr, 


dTr, 


dt 


(5) 


where the integral is over the mass of the H 2 -He envelope. 

Li n t is termed the “intrinsic” luminosity of the planet, 
that is, the net luminosity of cooling, equal to the total 
luminosity minus the “external” luminosity—the inter¬ 
cepted stellar flux—so that L int = L tot — L ext . Each 
of these components may be associated with a tempera¬ 
ture, most importantly, the effective temperature , defined 
by L int = 4,nR 2 aT^ ff . 

Energy is added to the envelope by radioactive decay, 
Lradio, and the cooling of the core, c v M core dT fj . re , which 
we model as an isothermal “hot rock,” where c v is the 
specific heat capacity. In this paper, unless otherwise 
specified, we define the core (or the composition of a 
purely “terrestrial” planet) as having an Earth-likc core- 
mantle structure of 32.5% Fe and 67.5% MgSiOa. 

The convective portion of the envelope is modeled 
with an eq uation of state for H 2 -He, as described in 
IHowe et al.l (l2014f) . with metals included as an ideal 
gas with the volume addition law. We include upper 
and lower boundary conditions to determine the struc¬ 
ture of the envelope. The lower boundary condition is 
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simply that the base of the envelope must be contin¬ 
uous in temperature and pressure with an isothermal 
core. For the upper boundary condition, we compute a 
grid of atmosphere models using the Coo lTLUSTY code 
dHubenv et alJl2003t iBurrows et al.ll200 6i , incorporating 
stellar irradiation, assuming complete redistribution and 
reradiation of the stellar flux over the surface. Interpo¬ 
lating in this grid provides a relation between effective 
temperature, entropy, and gravity. 

We decree that the base of the radiative atmosphere 
must be continuous in pressure and entropy with the 
isentropic envelope and use the entropy and gravity to 
determine the effective temperature at the radiative- 
convective boundary that we use to compute the cooling 
rate of the model. The radiative atmosphere also adds an 
additional amount, A R, to the radius of the planet com¬ 
puted by the structural code. This A R may be a signifi¬ 
cant fraction of the total radius of the planet throughout 
its lifetime, but especially at early times when the at¬ 
mosphere is “puffier”. In all cases, we report the total 
radius, R = R CO nv + A R. There is a discrepancy be¬ 
tween the total radius, which extends up to a pressure 
level of microbars and the observed transit radius, which 
is equivalent to cutting the atmosphere off at a lower al¬ 
titude corresponding to a pressure of tens of millibars 
(more specificall y a vertic al R ossel and optic al depth of 
r ~ 0.56: iLecavelier Des Etangs et al.ll2008l l. However, 
this is a secondary effect on the order of ~0.2 R® for a 
representative model in our set, and we find that it has 
a relatively small effect (< 10 %) on our fits of initial en¬ 
velope mass fractions to observed mass and radius data. 

For mass loss, we use an energy-limited scheme based 
on lNettelmann et al.l (120111) . which is approp riate for rel¬ 
atively lower irradition levels (j.Iin et al.ll2014fl and is com¬ 
mon in evolutionary modeling: 


= ^FxuyRxuv 

GM p K ti de 

Ktide = ( X - 2 £ + 2^) 
£ _ Rhui 
R xuv ’ 


( 6 ) 

(7) 

( 8 ) 


where Fxuv is the ionizing flux from the star, and Rxuv 
is the radius at which the atmosphere is optically thick 
to ionizing radiation, for which we use our total radius, 
R = Rconv + A R. Ktide is a correction for the tidal dis¬ 
tortion by the star of the planet (this is an approximate 
correction neglecting the dependence of tidal properties 
on the density profile of the planet), and e is an effi¬ 
ciency factor that we set by default to e = 0.1. Thus, we 
employ the approximation that a certain fraction of the 
XUV energy intercepted by the planet goes into strip¬ 
ping the atmosphere. This efficiency may vary under 
different circumstances and in parti cular may decrease 
at high levels of irradiation (I.Tin et alJl201 4). so we com¬ 
pute some models with a lower efficiency and keep this 
particularly in mind when trying to fit highly irradiated 
observed objects such as Kepler-lib and 11c (see Section 
ED- We consider the envelope to be lost if / env < 10 4 . 

For the XUV flux history, we use the model of 
iRibas et al.l (120051) for Sun-like stars, which gives a flux 


at 1 AU of 

Fxuv = 504 erg s ^ 1 cm~ 2 , t < 100 Myr (9) 

Fxuv = 29. Hg 1 ' 23 erg s ^ 1 cm~ 2 , t > 100 Myr ( 10 ) 

where tg is the age of the star in Gyr. iSanz-Forcada et ~ahl 
(l 2010 l) show that this XUV flux is approximately the 
same f or stellar types from F7 to M3. However, 
iShkolnik fc Barman! (120141) find a lower initial XUV flux 
for M dwarfs, so there is still some uncertainty in these 
models. 

We also consider the possibility that planets may form 
at relatively large distances where mass loss is negligi¬ 
ble and migrate in to close orbits where mass loss is 
high some time after formation, perhaps on the order 
of 100 Myr. Here, we are proposing migration by secular 
processes such as planet-planet scattering or the Kozai 
mechanism that may occur long after formation, followed 
by tidal circularization. We define a migration time , 
before which we set e = 0 (no mass loss) and after which 
we set e > 0 to simulate this late migration behavior in 
some models. By default, we set tM = 0, for the case 
where there is no migration from large distances, or such 
migration occurs during the planet formation period. 

The entropy and composition of the envelope fully de¬ 
termine the structure of the convective region, while for 
the radiative atmosphere, stellar irradiation must also be 
taken into account. In practice, to compute an evolution¬ 
ary track for a planet, we set an initial radius and fit an 
initial envelope entropy for the appropriate equation of 
state and atmosphere model. 

In principle, the initial radius should be the radius of 
the planet at the end of the accretion phase of planet 
formation. However, there are f ew es timates of “ini¬ 
tial radii” in the literature. iLee et all (I2014T ) modeled 
gas accretion onto solid cores and argued that the outer 
boundary of the planet during accretion should be the 
minimum of the Hill radius and the Bondi radius, on the 
order of ~40 R® for the planets in question. Ho wever, 
evolutionary models, e.g. iLopez fc Fortnevl (|2014H . more 
commonly set their initial conditions at a “hot start” of 
relatively high entropy or a “cold start” of relatively low 
entropy, resulting in initial radius of ~5-15 R®. While 
we set the initial radius directly instead of the entropy, 
we retain this strategy and use an initial radius range of 
5-20 R®, in part t o en sure speed and accuracy of compu¬ 
tation. In Section I7T7T1 we show that for R t > 10 R®, the 
precise initial radius matters very little to the subsequent 
evolution because the radius contracts very quickly. 

For evolution, we consider the envelope to be lost en¬ 
tirely when it reached f env < 0 . 01 % and halt the evolu¬ 
tion at that point. We also halt the evolution if the ef¬ 
fective temperature of the planet falls below 50 K, which 
is off our grid of atmosphere models, so a few models 
terminate before losing their entire envelopes. 

4. EXAMPLE INTERNAL STRUCTURE PROFILES 

Our code gives us information about the internal struc¬ 
tures of our planetary models, such as density, pressure, 
and temperature profiles of some sample models, and the 
evolution of those profiles over time. We plot these pro¬ 
files for two sample models in Figures[3] [H and[5l respec¬ 
tively. These models have parameters of initial M tot = 8 
M®, a = 0.30 AU, Ri = 10 R®, and Z/Zq = 30, and 
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Radius (R e ) 


Fig. 3.— Density profiles over time for a planetary model with 
initial parameters of Mtot = 8 M®, fenv = 0.01 and 0.05, Ri = 10 
R®, a = 0.30 AU, and otherwise using our default parameters. 



Fig. 4.— Pressure profiles over time for planetary models with 
initial parameters of Mtot = 8 M®, fenv = 0.01 and 0.05, Ri = 10 
R®, a = 0.30 AU, and otherwise using our default parameters. 

the specific models are those with initial f env = 0.01 and 
0.05. Both cases show negligible change in the rock-iron 
core over time, while the convective envelope gradually 
shrinks due to contraction in both cases (while also los¬ 
ing some mass), as indicated by the increasing density 
near the core-envelope boundary. The temperature of 
the core and the temperature profile of the envelope also 
decrease significantly over time. 

5. EXAMPLE MASS LOSS RATES 



Fig. 5. — Temperature profiles over time for planetary models 
with initial parameters of Mtot = 8 M®, fenv = 0.01 and 0.05, 
Ri = 10 R®, a = 0.30 AU, and otherwise using our default param¬ 
eters. 



Fig. 6. — Plot of the mass loss rates of mini-Neptune models 
with a = 0.30 AU, and our “default” parameters of Ri = 10 R®, 
Z/Z® = 30, and e = 0.10. Orbital distances of 0.07 AU, 0.10 AU, 
0.14 AU, 0.20 AU, and 0.30 AU are plotted with an initial envelope 
fraction of fenv = 0.05 and initial total masses of 2, 4, and 8 M®. 

In Figure [51 we provide mass loss rates over time for 
some sample planet models. We compute models with 
initial total masses of 2 M®, 4 M®, and 8 M®, and orbital 
distances of 0.07 AU, 0.10 AU, 0.14 AU, 0.20 AU, and 
0.30 AU, with an initial envelope fraction of f en v = 0.05 
and with initial radii of 10 R®. Given the relation M oc 
FX K t v d R •> we ex P ec t the mass loss rates to fall off roughly 
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in line with the XUV flux, with smaller influences from 
the decreasing radii and tidal corrections. This is indeed 
what we see, with the mass loss rates being flat or varying 
relatively little at t < 100 Myr, the time during which 
XUV flux is constant, decreasing for models that undergo 
significant contraction during this period. At t > 100 
Myr, the mass loss rates fall off at on the order of the 
same rate as the XUV flux. We also see mass loss rates 
decreasing with increasing orbital distances due to the 
lower XUV flux. 



Fig. 7.— Plot of mass loss rate versus orbital distance at an 
age of 2 Gyr for mini-Neptune models. We include models with 
Ri = 10 R, 0 , Mi = 8 M®, and f env ,i = 0.01, 0.02, 0.05, and 0.10. 

As a further illustration, in Figure [7] we plot the mass 
loss rate versus distance at 2 Gyr for select models, where 
we set Mi = 8 M®, and we plot f env ,i of 0.01, 0.02, 0.05, 
and 0.10. Here the behavior is as expected with M oc 
a~ 2 /K t ide from about 0.10 AU to 0.30 AU. The turnoff at 
the closest distances, 0.07-0.10 AU, occurs because most 
of the envelope has been lost at this point, decreasing 
the radius of the planet and, thus, the mass loss rate. 

6 . PARAMETER STUDY 

In this section, we study in depth the effects of six 
parameters on the models produced by our code: initial 
total mass (Mi), initial envelope fraction (f en v,i)i orbital 
distance (a), initial radius (Ri), metallicity (Z/Zq), and 
mass loss prescription (e and £m)- Figures 151(171 show the 
results of our calculations. 

6 .1. Dependence on Initial Mass 

The mass of the solid (terrestrial-composition) core 
comprises most of the total mass in our models. In 
Figures |U and [9] we plot models with M* = 2, 4, and 
8 M®. In Figures [TO] [Til [12] and [5] we plot theo¬ 
retical mass-radius relations for our models at differ¬ 
ent times over an initial mass range of 0.5-10 M®. 
We plot these models with a mass-radius relation for 


a purely terrestrial composition ( black line) and obser ¬ 
vati onal data (b l ack p oints) from iDressing et al.1 (|2015fl 
and I Weiss et all (|2015fl . Note that our theoretical mass- 
radius curves are at fixed orbital distances, which are not 
the same as those of the observed planets. The observed 
planets have different ages and orbital distances than our 
model curves and serve simply to illustrate the range of 
planets we can model by comparison. 

For comparison, low-density exoplanets have been ob¬ 
served with masses as small as ~2 M®, while the largest 
terrestrial planets are ~8 M® (see Fig. [I]), although this 
limit may be significantly greater for icy planets. Also, it 
is generally thought that 10 M® is the critical core mass 
for runaway gas accretion, although there is so me un¬ 
certainty in thi s, as sugg ested by iRafikovI (j2006T ). More 
recently. iLee et all (|2014H . found that 10-M® rocky cores 
readily undergo runaway gas accretion regardless of or¬ 
bital distance, but they speculate as to mechanisms that 
might prevent this and lead to the observed preponder¬ 
ance of super-Earths. Therefore, 2-10 M® may be con¬ 
sidered a plausible range of core masses. 

We find that planets with lower initial masses experi¬ 
ence more mass loss due to lower gravity, although the 
amount of mass loss is more sensitive to orbital distance 
(see Section R751) . The lowest-mass planets lose their en¬ 
velopes completely, although the limit of this effect also 
varies with orbital distance. However, for higher masses, 
the mass-radius curves are remarkably insensitive to ini¬ 
tial mass. We also find that most of our models have “fi¬ 
nal” radii of 2-4 R®, the range where the greatest number 
of Kepler planet candidates have been observed. 

6 .2. Dependence on Envelope Mass Fraction 

In Panel (a) of Figures [8] and [9] and all of Figure [101 we 
see the effect of varying the initial envelope mass fraction, 
fenv,i on our model evolution, with f en v,i = 0.01, 0.02, 
0.05, and 0.10. Notably, for these and all of our models, 
we find that the radius evolution nearly ceases at Gyr 
ages because both cooling and mass loss are much slower 
then than at earlier times. This removes one important 
degeneracy when trying to fit theory to observed planets: 
for most observed planets, which are usually older than 
1 Gyr, their properties will have very little dependence 
on their age. 

Meanwhile the “final radii” (or radii at Gyr ages) can 
vary significantly depending on the masses of both the 
solid core (which constitutes most of the total mass) and 
the gaseous component at formation. For the most part, 
this applies to f en v,i in particular; the mass-radius curves 
are remarkably flat above a certain mass, with variations 
in radius of < 0.1 R® over a wide range of masses. In any 
case, the difference in radius at Gyr ages can be multiple 
Earth radii with different initial envelope fractions. 

We also find that the amount of mass loss as a fraction 
of the initial envelope mass varies relatively little with 
f env , especially at higher masses (e.g., ~8 M®), where 
it is effectively negligible. However, this factor is more 
important at lower masses. 

6.3. Dependence on Orbital Distance 

In Panel (b) of Figures [5] and GO and Panel (a) of Fig¬ 
ures HH we present models with a range of orbital dis¬ 
tances: 0.07, 0.10, 0.14, 0.20, and 0.30 AU. This is the 




















Parameter Study of Radius v. Time for Mini-Neptune Models 
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Fig. 8 .— Parameter study of radius evolution of mini-Neptune models. In all cases, models with initial masses of 2, 4, and 8 M® are 
shown. Panel (a): study of initial envelope fraction with f e nv,i = 0.01, 0.02, 0.05, and 0.10. Panel (b): study of orbital distance with a = 
0.07, 0.10, 0.14, 0.20, and 0.30 AU. Panel (c): study of initial radius with Ri = 5, 10, 15, and 20 R®. Panel (d): study of metallicity with 
Z/Z® = 1, 3, 10, and 30. 


distance range (stellar flux range) in which most of the 
mini-Neptunes with significant envelopes have been dis¬ 
covered. From these, we find that the shape of our theo¬ 
retical mass-radius curves (at least at Gyr ages) is similar 
across our models, but shifted left or right, to lower or 
higher masses, depending upon orbital distance. This is 
most readily seen in Panel (a) of Figure |TT] 

For a fixed f env ,i , the mass-radius curve at 2 Gyr is ap¬ 
proximately flat at high masses, while at lower masses, it 
reaches a “turnoff point” that depends on the orbital dis¬ 
tance, at which the radius rapidly decreases over a range 


of roughly 2 until it reaches the terrestrial core (com¬ 
plete envelope evaporation). For closer orbital distances, 
this turnoff point occurs at higher masses, approximately 
1-2 higher for each doubling in stellar irradiation. (It 
also occurs at somewhat lower masses for higher f en v,i-) 
The radius at 2 Gyr is remarkably consistent over a 
range of masses, and, therefore, over a range of enve¬ 
lope masses and surface gravities as well. To investigate 
how this relates to the amount of mass loss, in Panel (b) 
of Figure [9] we plot the envelope mass evolution, that 
is the fraction of the mass of the envelope remaining, 
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Parameter Study of Mass v. Time for Mini-Neptune Models 



Fig. 9.— Parameter study of envelope evaporation of mini-Neptune models. In all cases, models with initial masses of 2, 4, and 8 M® 
are shown. Panel (a): study of initial envelope fraction with f e nv,i = 0.01, 0.02, 0.05, and 0.10. Panel (b): study of orbital distance with 
a = 0.07, 0.10, 0.14, 0.20, and 0.30 AU. Panel (c): study of initial radius with Ri = 5, 10, 15, and 20 R®. Panel (d): study of metallicity 
with Z/Zq = 1, 3, 10, and 30. 


versus time for selected models (all with f env ,i = 0.05) 
at a range of orbital distances. Many of the resultant 
curves represent models that converge to approximately 
the same radius at Gyr ages, but nonetheless undergo 
very different amounts of mass loss, losing anywhere from 
5% to nearly 100% of their envelopes. Thus, the collapse 
of the envelope due to cooling is significantly more de¬ 
pendent on orbital distance than on other parameters, 
regardless of the density and pressure in the envelope. 
The pattern of mass loss increasing with decreasing or¬ 
bital distance is as expected from our model, where the 


mass loss rate scales as a _2 /A' t i C j e (where K t ide is the 
tidal correction). 

In Figure [T2] we extend our plot of mass-radius curves 
to 20 M® and orbital distances of 0.02 AU to 0.30 AU. 
Here, we find that the above trend in the shape of the 
mass-radius curve extends to both of these parameters, 
to some extent. On the high-mass end, the mass-radius 
curves remain approximately flat to 20 M®. Mean¬ 
while, the previously-produced behavior of a turnoff 
point (where, for lower masses, the radius decreases 
rapidly) continues to an orbital distance of 0.05 AU, 
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Fig. 10.— Study of mass-radius curves with different f e nv,i at 10 Myr and 2 Gyr for mini-Neptune models. In all cases, models with 
a = 0.30 AU, and Ri = 10 R®, Z/Zq = 30, and initial total masses of 0.5-10 M® are plotted. Panel (a): f en v,i = 0.01. Panel (b): 
fenv,i = 0.02. Panel (c): / en t),i = 0.05. Panel (d): f enV) i = 0.10. 


where it occurs at a higher mass and leaves higher-mass 
bare cores. However, for 0.03 AU, complete evaporation 
of the envelope occurs almost all the way up to 20 M®. 
This is consistent with the fact that no planets with likely 
gaseous envelopes are observed with stellar fluxes greater 
than ~700 F®, equivalent to about a = 0.04 AU around 
a sun-like star. 

To further illustrate the effects of orbital distance, in 
Figures [13] and [Ufl we plot several quantities for represen¬ 
tative mini-Neptune models versus orbital distance over 
a range of 0.07-0.30 AU. In each case, we plot the states 
of the models at an age of 2 Gyr, and we plot models with 
fenv.i = 0.05 and M» = 2, 4, 8, and 16 M®, otherwise 
using our default parameters. 

Figure [T3] shows a plot of radius at 2 Gyr versus orbital 
distance. We again see that the radius is relatively flat 
with orbital distance at lower levels of irradiation, but 
rapidly falls off at closer orbital distances because the 
envelope has mostly or entirely evaporated. This turnoff 
occurs at 0.10 AU for the 16-M® models, but farther out 
for lower masses to the point where it occurs outward of 
0.30 AU for the 2-M® models. 

Figure [14] plots the fraction of the envelope mass re¬ 
maining after 2 Gyr versus orbital distance. This quan¬ 
tity follows a general pattern of increasing from zero or 


near zero to a high level over a small factor of 2-3 in or¬ 
bital distance. For example, the envelope retention fac¬ 
tor for the 16-M® models increases from zero at 0.03 AU 
to ~0.85 at 0.10 AU. At the same time, this same shift 
occurs much farther out for the 2-M® models, with an 
envelope retention factor of ~0.03 at 0.20 AU and ~0.22 
at 0.30 AU, retaining a majority of the initial envelope 
only at larger distances. 

6.4. Dependence on Initial Radius 

The initial radius of a given planet at the end of accre¬ 
tion is unclear. Therefore, we explore models of planets 
with H 2 -He envelopes with initial radii in the range from 
22.6 iuE near the pr obab le maximum size o f exo plan- 
ets of any mass dMordasini et al.ll2012at lAnderson et al.l 
[MU . down to a radius of ~2-3 R©0 

In Panel (c) of Figures [5] and O and Panel (b) of Figures 
El we present models with a range of initial radii: 5, 10, 
15, and 20 R®. Ri = 10 R® produces the widest range 
of well-defined models, which is why we have used it as 
a default throughout our parameter study. 

2 Equivalent to 2 Jupiter radii. 

3 The actual limit is T e ff = 50 K, the boundary of our atmo¬ 
sphere grid, which corresponds to a radius of 2-3 R® for some of 
our representative models. 
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Fig. 11.— Parameter study of mass-radius curves at 10 Myr 
and 2 Gyr for mini-Neptune models. In all cases, models with 
fenv,i = 0.05 and initial total masses of 0.5-10 Mj are plotted. 
Panel (a): study of orbital distance with a = 0.07, 0.10, 0.14, 0.20, 
and 0.30 AU. Panel (b): study of initial radius with R, = 5, 10, 
15, and 20 R®. Panel (c): study of metallicity with Z/Zg = 1, 3, 
10, and 30. 


By the assumptions of our model, mass loss is pro¬ 
portional to the amount of starlight intercepted by the 
planet and is inversely proportional to the gravitational 
potential at the XUV photosphere, thus making it pro¬ 
portional to R?/Ktide- Therefore a larger initial radius 
leads to much greater mass loss and could even lead to 
complete evaporation of the envelope. ( 1/Kude also in¬ 
creases with radius.) 

Notably, however, we find that the initial radius of the 
models makes only a small difference in either radius or 
mass loss at Gyr ages. This is because of rapid radius 
contraction at early times, causing models differing in 
initial radii to converge or nearly converge by an age of 
100 Myr, and often much earlier. This early convergence 
in radius also results in very small differences in mass 
loss, which is primarily dependent on radius for equal- 
mass planets. The exceptions to this rule are models 
with low initial total mass, which fall below the turnoff 
point. This removes one source of degeneracy in trying 
to fit observational data. 


6.5. Dependence on Metallicity 


Mass-Radius Relations at Different Orbital Distances 



Fig. 12.— Mass-radius curves for planetary models with large 
masses of 10-50 M®, initial fenv = 0.01, 0.02, 0.05, 0.10, and 
0.20, initial radii of 10 R. 0 , and orb ital distances o f 0.30 AU, 
plotted with obse rvational data from [Dressing et al. (2015[) and 
Weiss ct al. (2015). Curves at 100 Myr and 5 Gyr are shown. 



a (AU) 


Fig. 13.— Plot of radius versus orbital distance at an age of 2 
Gyr for mini-Neptune models. We include models with Ri = 10 
R -0 5 fenv,i — 0.05 M 0 , and Mi = 2, 4, 8 , and 16 M 0 . 

In Panel (d) of Figures[8]and[9l and Panel (c) of Figures 
EU we present models with a range of envelope metallici- 
ties: lx, 3x, lOx, and 30x solar metallicity. For compari¬ 
son, atmospheres of super-Earths and mini-Neptunes are 
expected to be significantly enriched compared with so¬ 
lar. Specifically, Uranu s and Nep tu ne ar e enriched to 
~50 times solar values dHelled et al.l 1201 ill . While it is 
not obvious a priori that close-in planets should be sim- 
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Fig. 14.— Plot of envelope retention factor versus orbital distance 
at an age of 2 Gyr for mini-Neptune models. We include models 
with Ri = 10 R®, fenv : i = 0.05, and Mi = 2, 4, 8 , and 16 M®. 

ilarly enriched, it is a reasonable assumption since ob¬ 
servations of transiting super-Earths and sub-Neptunes 
support either high-altitude cl ouds and/or a high me an 
molecular weight atmosphere (|Howe fc Burrowsll2012f) . 

Varying the metallicity of the models has some effect, 
though not as large as that of initial total mass or orbital 
distance. A higher metallicity increases the mean molec¬ 
ular weight of the atmosphere, which should, in principle, 
increase the density. However, we find in our simulations 
that this is a small effect on the convective region and 
that the greater opacity, which results in a thicker radia¬ 
tive layer and slower contraction time, results in a larger 
overall radius at a given age for higher metallicity (de¬ 
spite higher mass loss). Thus, a higher metallicity results 
in a larger radius at Gyr ages. The difference in final ra¬ 
dius between lx and 30x solar metallicity is about 0.40 
R®, although this effect is more dramatic at early times. 

Because of this smaller variation, and because the at¬ 
mospheres of mini-Neptune-type planets are expected to 
be highly enriched in general (>30x solar), metallicity is 
not a critical parameter to fit observational data. 

The effect of metallicity on mass loss is similar to the 
effect of initial envelope fraction. It is relatively unimpor¬ 
tant compared with other parameters, and it is especially 
negligible for the higher mass objects in our models with 
Mi = 8 M®. 

6 .6. Dependence on Mass Loss Prescription 

While the relationships between mass loss rates and 
physical parameters are straightforward, the processes 
that lead to mass loss in gase ous exoplanets are uncer¬ 
tain (iMurrav-Clav et aD 12009 1. However, especially for 
a specific model, the mass loss behavior may be parame¬ 
terized by the efficiency—the fraction of intercepted XUV 
energy that goes into mass loss. In most of our models 
we set the efficiency at e = 0.10. We also have the option 
of having mass loss initially switched off and then turned 



Fig. 15.— Comparison of different mass loss prescriptions in plots 
of radius versus time for mini-Neptune models. Plotted are e = 0 
(no mass loss), e = 0.05, e = 0.10, and e = 0.10 with the mass loss 
switched off until 100 Myr. We include models with Ri = 10 R, 0 , 
a = 0.30 AU, Mi = 4 M®, and f e nv,i — 0.01, 0.02, 0.05, and 0.10. 



Time (Myr) 


Fig. 16.— Comparison of different mass loss prescriptions in 
terms of the mass fraction of the original envelope remaining versus 
time for mini-Neptune models. Plotted are e = 0 (no mass loss), 
e = 0.05, e = 0.10, and e = 0.10 with the mass loss switched off 
until 100 Myr. We include models with Ri = 10 R, 0 , a = 0.30 AU, 
Mi = 4 M 0 , and f e nv,i = 0.01, 0.02, 0.05, and 0.10. 

on at some point in time, tM■ This simulates a planet 
that forms at a relatively large distance, where mass loss 
is negligible, and then migrates to a close-in orbit where 
it undergoes higher mass loss at time tM- 
In Figures [15] and [lGl we plot radius and mass, respec¬ 
tively, versus time for planet models with four different 
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mass loss prescriptions: our default e = 0.10, e = 0.05, 
e = 0 (no mass loss), and a fourth model with e = 0.10 
switched on at tM = 100 Myr. This type of prescription 
can bracket the actual, smoother variation in mass loss 
rates for migrating planets. We plot models with f env = 
0.01, 0 . 02 , and 0.05 and a = 0.30 AU, and we otherwise 
use our default parameters. 

With no mass loss, the only influence on the radius is 
cooling, and the final radii are significantly larger than 
with mass loss, by up to about an Earth radius. No¬ 
tably, without mass loss, the envelopes of low-mass plan¬ 
ets never shrink to small radii, and so they are not consis¬ 
tent with low-mass, low-density planets such as Kepler- 
11b, -11c, and -Ilf. Some mass loss is required to fit 
our models to observations in this part of the parameter 
space. 

Setting e = 0.05 leads to half as much mass loss as 
our default models and makes up approximately half the 
difference in radius. The models in which the mass loss 
switches on only after 100 Myr also experience about half 
as much mass loss, both in terms of radius and actual 
mass. Late migration due to planet-planet scattering or 
the Kozai mechanism beginning at t > 200 Myr would 
have a much greater impact on mass loss because the 
XUV radition falls off as t~ 123 after 100 Myr. 

In Figure 1171 we plot mass-radius curves over a range 
of 0.5-10 M® for the same four cases. For this figure, 
we plot the mass-radius curves at 10 Myr and 2 Gyr. 
We also plot a range of orbital distances: 0.07, 0.10, 
0.14, 0.20, and 0.30 AU (colored lines). We find that 
for all of our mass loss prescriptions except the case of 
no mass loss, several of the 2 Gyr curves pass through 
the low-density (low mass and high radius) planets on 
the plots, but they do not necessarily correspond to the 
same orbital distances. We investigate fitting our models 
to specific observed planets in Section 0 

7. APPLICATION TO THE KEPLER-11 SYSTEM 

Kepler-11 is perhaps the best- studied planetary sys¬ 
tem in terms of the structures of its planets. It is 
a densely-packed system with six known planets with 
se mi-maj or ax es ra nging from 0.091 AU to 0.466 AU 
dLissauer et alJ[ 20 TH) . The masses of five of these planets 
have been precisely measu red with transit timing varia¬ 
tions dLissauer et_a l. 20Uf) and , more recently, with ra¬ 
dial velocities (IWeiss et al.112015 ). The measured orbital 
distances, masses, and radii for these planets are shown 
i n Table [D 

iLopez et al.l (|2012f) studied the compositions of these 
planets using an evolutionary model with a H 2 -He en¬ 
velope surrounding an Earth-like core. However, these 
models were fit to the previous masses measured from 
TTVs, and the improved RV measurements invite a 
reanalysis of these planets. The measured masses for 
Kepler-lld, -lie, and -Ilf changed by a relatively small 
amount between the TTV and RV observations, but the 
masses (and densities) of Kepler-llb and -11c were re¬ 
vised down dramatically, by more than a factor of 2 (with 
significantly smaller error estimates), so their estimated 
compositions will be very different. 

We can change the stellar type in our model, which 
we have done for the planets in the Kepler-11 system. 
Changing the central star in our model has two effects. 
First, the total thermal irradiation changes according to 


a~ 2 R 2 Tf. Secon d, different st ars will have different ul¬ 
traviolet fluxes. iLouez et al .1 ( 120121 ) assumed that all 
planet host stars have the same XUV flux evolution. 
Ho wever, ot her models _suggest different prescriptions. 
iShkolnik fc Barmanl ( 2014 ) find that the power-law de¬ 
cay of XUV flux only holds at > 1 Gyr for M dwarfs, and 
these stars begin with a lower XUV flux, which remains 
constant longer. However, Kepler-11 is a G dwarf, so this 
is not a concern. 

A lower irradiation level would allow for more rapid 
cooling and contraction, leading to less mass loss, while 
a lower XUV flux also results in less mass loss. There¬ 
fore, around smaller stars, low-mass planets could retain 
gaseous envelopes longer and at closer distances than 
around larger stars. 

We find in general that the most important parame¬ 
ters to fit the data are initial mass (Mi), initial envelope 
fraction ( f env ,i ), and orbital distance (a, which is known 
from observations), while we can choose Ri = 10 R®, 
Z/Z Q = 30, and t > 2 Gyr without loss of generality, 
although we use a larger initial radius for the lower mass 
planets, Kepler-llb, -11c, and -Ilf, which produces bet¬ 
ter convergence with our codeQ 

Kepler-lld and -lie are the easiest objects to inter¬ 
pret, having high enough masses and large enough or¬ 
bital distances to retain significant amounts of hydrogen 
and helium over several Gyr. Given how closely the ini¬ 
tial envelope fraction correlates with radius at Gyr ages 
in our model set, it is possible to determine f env ,i , or 
the initial envelope mass, M env ^ fairly precisely. The 
uncertainty in the observed mass of 5%-10% introduces 
an uncertainty of similar size into M env ^, but the fit of 
M enVi i is positively correlated with the assumed mass of 
the model, so this uncertainty is diminished in f env ,i- 
Also, the roughly 2% uncertainty in radius introduces a 
5%-10% uncertainty into both M enV} i and f en v,i- 

For Kepler-lld, we find an initial H 2 -He content of 
0.25 ± 0.03 M®, corresponding to f env ,i = 3.7 ± 0.4%. 
Mass loss in this model is modest. At an age of 4 
Gyr, which will be very similar to the result at the 
system’s age of 8.5 Gyr, the resultant envelope frac¬ 
tion is fpm.ii = 2 .6 ± 0.3%. This is in contrast with 
ILopez et al.1 (120121) . who find a present-day envelope frac¬ 
tion of fenv = 8.4+24% f° r Kepler-lld (with older, but 
similar mass measurements). 

For Kepler-lle, we find an initial H 2 -He content of 
0.68 ± 0.06 M®, corresponding to f en v,i = 8.3 ± 0.4%. 
At an age of 4 Gyr, which will be very similar to the 
result at the system’s age of 8.5 Gy r, the resultant en- 
velope fraction is f env = 6.5 ± 0.5%. ILopez et al.1 (120121 ) 
find a present day f env for Kepler-lle of 17.2142% for 
Kepler-lle. For Kepler-lld and -lie, they estimate the 
initial envelope fraction, f en v,i at about 28%. Both of our 
estimates for present-day envelope fractions are lower by 
a factor of three. 

Part of this discrepancy could conceivably result from 
our treatment of the radiative atmosphere. (Specifically, 
the thickness of the atmosphere, since the opacities affect 
the cooling rate and are a separate problem.) We deter¬ 
mine an upper limit to the magnitude of this effect by 

4 Note that the evolutionary models used bv ILopez et al.l 112012! ) 
employed Ri = 10 — 15 R®, except for Kepler-llb, for which R, ~ 
30 R® was used. 
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Mass-Radius Relations Under Different Mass Loss Prescriptions 



Fig. 17.— Mass-radius curves at 10 Myr and 2 Gyr for mini-Neptune models with different mass loss prescription s, pl otted with a 
terrest rial, Earth-like composition model (black line), and observational data (black points) from fDressing et al.l <120151) and I Weiss et all 
(120151) . In all cases, models with Ri = 10 R®, f e nv,i = 0.05, orbital distances of 0.07, 0.10, 0.14, 0.20, and 0.30 AU, and initial masses of 
0.5-10 are plotted. Panel (a): e = 0 (no mass loss). Panel (b): e = 0.05. Panel (c): e = 0.10. Panel (d): e = 0.10 with the mass loss 
switched off until 100 Myr. 


TABLE 1 

Mass and Radius Observations of the Kepler-11 Planets 


Planet 

Orbital Distance (AU) 

Radius (R®) from 
Lissauer et al. (2013) 

Our adopted mass (M®) from 
Weiss et al. (2015) 

Mass (M®) adopted by 
Looe^et al. (2012) 

Kepler-lib 

0.091 

1.80 ±0.04 

1 07+0-43 

^-2.0 

Kepler-11c 

0.107 

2.87 ±0.05 

9 ftq+0.80 

z.uo_063 

13.511? 

Kepler-lld 

0.155 

3.12 ±0.06 

IJ q^J + 0.03 
f .OO-i 16 

e-iti! 

Kepler-lie 

0.195 

4.19 ±0.08 

o i q+ 0.77 
°- iy -0.66 

8-4±?;g 

Kepler-Ilf 
Kepler-llg 

0.250 

0.466 

2.49 ±0.05 

3.33 ±0.07 

1 Q4+0- 32 

0.88 

5 20+27.6 

9 o-|-2.2 

1.2 

N/A 


fitting a model to Kepler-lld with no A R at all—that is, 
no radiative atmosphere. With this fit, we find an initial 
envelope fraction of f e nv,i = 5.4 ± 0.3% and a present 
day envelope fraction of f env ,i = 4.6 ± 0.2%, st ill a fac¬ 
tor of two lower than the results of lLonez et al.l 1 2012 1. 
Therefore, our treatment of the radiative atmosphere is 
likely not a large source of error in our analysis. 

It is also possible that the difference in computed enve¬ 
lope fractions is partially due to the choice of equation of 
state, specifically, the use of the volume addition law for 


the envelope, and the density of the solid core (although 
the latter effect is likely small), both of which could affect 
the determination of the relative core-envelope fractions. 
However, the most important factor is likely the cooling 
rate of the models, based upon both the radiative prop¬ 
erties of our atmospheric boundary conditions and the 
specific heat of the core, which contribute to the overall 
cooling at the same order of magnitude. Faster cooling 
results in more contraction and a larger envelope fraction 
to achieve the observed radii. Further work is needed to 
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determine the magnitude of these effects and the most 
correct boundary condition. 

The other three planets, Kepler-llb, -11c, and -Ilf, 
are much more difficult to interpret. Kepler-llb and 
-11c were both revised down in mass by more than a 
facto r of two compared with the measurements used by 
iLopez et al.l (120121) . who already find a high mass loss 
rate, requiring a large initial envelope fraction. A lower 
mass, and thus a lower gravity, will increase this problem 
further. 

In the case of Kepler-llb, the current radius of 1.80 R® 
is too small to explain without evaporation of most, but 
not all, of its envelope. With such a low mass (1.87 M®) 
and a relatively high level of irradiation (a = 0.091 AU), 
the atmosphere puffs up enough that it is difficult to find 
a solution that achieves the observed radius. Further¬ 
more, those solutions that do reach the observed radius 
are all close to evaporating completely by the present 
day, even if the evaporation efficiency is set an order of 
magnitude lower than expected, to e = 0.01, or if evapo¬ 
ration starts late due to migration ( Im = 300 Myr) due to 
planet-planet scattering or the Kozai mechanism. Thus, 
at expected evaporation rates, most models are likely 
to lose their envelopes completely, which is inconsistent 
with observations. Models with much higher starting 
masses (tens of M®) stay larger than the observed ra¬ 
dius. 

Because of its high irradiation (a = 0.107 AU), Kepler- 
11c is likely to be subject to the same problems as Kepler- 
llb in finding a solution that retains its envelope to the 
present time, although its mass is higher, which aids 
in envelope retention. Because its radius is larger, the 
“water-world” solution does not hold; its low density re¬ 
quires some hydrogen. For this object, we still find that 
all of our models lose their envelopes completely for our 
standard mass loss rate with an efficiency of e = 0.10. 
However, we find a solution of f en v,i ~ 1.5% when the 
mass loss efficiency is set to e = 0.01. For this solution, 
the present-day envelope mass fraction is f env ~ 0.25%, 
indicating a large amount of mass loss. 

However, this seems to be an improbably low mass loss 
efficiency, so there may be a different mechanism at work. 
If we leave e = 0.10, but set tM = 300 Myr, indicating 
that the planet moved to its current orbit from farther 
out at an age of 300 Myr, we find a similar solution with a 
somewhat larger envelope of f enVt i ~ 0.025 and a present- 
day envelope fraction of f env ~ 0.7%. We fail to produce 
such a result with a migration time significantly earlier 
than 300 Myr. 

In light of this, we construct an alternative solution 
that explains the observations of both Kepler-llb and 
Kepler-llc. In this scenario, Kepler-llb is a “water 
world” that formed beyond the snow line with a compo¬ 
sition consisting mostly of water, including an ice core, 
and either no hydrogen or helium or an envelope that 
has since evaporated completely, leaving a thick steam 
atmosphere. (T his model is an estima te b ased on our 
previ ous work in iHowe fc~Bu rrows 2 0121 and iHowe et al.1 
2014.) Kepler-llc also formed below the snow line as 
a mini-Neptune with an initial envelope mass fraction 
of ~2.5%. The two planets then moved to their cur¬ 
rent orbits some time after they formed, on the order 
of 300 Myr. In their current orbits, Kepler-llb lost any 
H 2 -He envelope that it had, leaving its thick steam atmo¬ 


sphere, and Kepler-llc lost most of its envelope, reach¬ 
ing its current composition of f env ~ 0.7%. The water 
world explanation is the only model we find for Kepler- 
llb that fits the observations within observational un¬ 
certainty. However this fit is provisional, and additional 
research is needed to determine its plausibility. 

For Kepler-llf, for the standard mass loss prescrip¬ 
tion of e = 0.10, we find an initial envelope fraction at 
the high-mass end of the uncertainty range of f en v,i = 
2.7 ± 0.2%. The models that best fit the observations 
fall near the edge of our grid of atmospheres, but we can 
extrapolate to a significantly higher fraction f en v,i ~ 4% 
near the middle of the uncertainty range. The present- 
day envelope fraction is on the order of f env ~ 0.25%, 
indicating significant mass loss. Both the initial and fi¬ 
nal envelope mass fractions in our analysis are an order 
of ma gnitude smaller than those found by iLonez et al.1 
llOll). 

The latest available data also include a mass mea¬ 
surement for Kepler-llg, which is 5.20 M®. The 2cr 
upper mass limit for this object is 27.6 M®. Assum¬ 
ing a true mass near this preliminary measurement, we 
find an initial envelope mass fraction for Kepler-llg of 
f env ,i = 2.7±0.3%. Because Kepler-llg is at a relatively 
large distance of 0.466 AU, the mass loss is small, with 
a present-day envelope fraction of f env = 2.6 ± 0.3%. 

8. EXAMPLE EVOLUTIONARY PROGRESSIONS OF 
TRANSIT SPECTRA 

Given the masses, radii, temperatures, and atmo¬ 
spheric compositions of our models over time, we can 
compute th e evolution of thei r trans it spectra using the 
methods of IHowe fe Burrowsl (120121 1. In Figure fT8l we 
present the evolution of two sample spectra. These 
models have parameters of initial M tot = 8 M®, initial 
f env = 0.05, Ri = 10 R®, and Z/Zq = 10, and the spe¬ 
cific models are those with a = 0.10 (left panel) and 0.30 
(right panel). 

For the purpose of these spectral models, we consider 
an isothermal radiative atmosphere. We also model a 
cloudless atmosphere because that is most consistent 
with our evolutionary models. Observations of atmo¬ 
spheres of mini-Neptune-type p lanets remain ambiguou s 
and may support a haze layer (jHowe fe Burrowsll2012il . 
as has been observed for giant planets. However, we note 
that our results were relatively little affected by metallic- 
ity, so the effect of clouds or hazes on evolution is likely 
to be similarly small. Moreover, this cloudless model is 
similar t o that used in oth e r evol utionary models such as 
those of lLopez fe Fortnevl (|2013l) . 

The temperature of this isothermal model is equal 
to that of our full atmosphere model at a pressure 
level corresponding with a vertical Rosseland opti¬ 
cal depth of r = 0.56, which is approximately the 
mean optical de pth pro bed by the transit spectrum 
(jLecavelier Des Etangs et al.l 120081) . Real atmospheric 
temperatures change little over a Gyr, primarily because 
of stellar irradiation, so the main sources of change in the 
spectra are the decrease in the overall transit depth due 
to contraction and the increase in gravity. As a result, 
the evolution of the spectra is mostly a flattening, and 
not a change, in the features. 

However, between the two orbital distances, we see 
significant differences. First, the atmosphere tempera- 



























16 


C 

0 ) 

o 

t-, 

<u 

cu 


si 

Oh 

(U 

Q 


• pH 

w 
c 
(0 
l H 


0.8 

0.7 

0.6 

0.5 

0.4 

0.3 

0.2 

0.1 



Wavelength (jum) 


Wavelength (/um) 


Fig. 18.— Transit spectra over time for a planetary model with initial parameters of Mtot = 8 M®, fenv = 0.05, and Ri = 5 R®, 
and Z/Z® = 10. We assume a cloudless atmosphere (predicted by the simplest theoretical treatment of the given modeled atmosphere 
temperatures of 500-900 K for all models). Left panel: a = 0.10 AU; right panel: a = 0.30 AU. 


tures provided by our code as described above are very 
different—nearly 800 K at 0.10 AU and only about 450 K 
at 0.30 AU. The higher temperature results in a higher 
scale height and a deeper atmosphere. Also, at the higher 
temperature, we see prominent Na and K alkali metal 
lines in the spectrum that we do not see at lower tem¬ 
peratures, and we see strong CO features and weak CH 4 
features at 0.10 AU, whereas the opposite is true at 0.30 
AU. 

9. CONCLUSIONS 

We have constructed exoplanet evolution models for 
super-Earths and mini-Neptunes incorporating thermal 
cooling and XUV-driven mass loss over a wide range 
of parameters. With these models, we have conducted 
a parameter study of the effects of initial conditions: 
mass, radius, envelope mass fraction, orbital distance, 
and metallicity, as well as the particular prescription for 
mass loss on super-Earth and mini-Neptune evolution. 

We find that the most important parameters to pro¬ 
duce the observed variation in mini-Neptune type exo¬ 
planets are orbital distance, initial envelope mass frac¬ 
tion, and initial total mass, while age, metallicity, and 
initial radius make relatively little difference in the 
present-day observable properties of many such exoplan¬ 
ets, initial radius being the least important. Our gen¬ 
eral results show that a large subset of our models have 
present-day radii of 2-4 R®, consistent with the radius 
distribution of Kep ler pla net ca ndidates, which exhibits 
a peak at 2-3 R® (iBurke et al.ll2014f) . These radii are 
reached within ~2 Gyr, with only slow radius evolution 
at later times. 


Moreover, we find that the present-day radius of a 
planet is a very good proxy for its initial envelope mass 
fraction above a certain “turn-off mass,” below which 
mass loss dominates, and planets rapidly become “ter¬ 
restrial.” This turn-off mass is dependent upon orbital 
distance and ranges from about 10 M® at 0.05 AU to 
about 2 M® at 0.30 AU, shifting by approximately 1-2 
M® for each factor of 2 in stellar irradiation. At higher 
masses, up to 20 M®, our mass-radius curves are ap¬ 
proximately flat and mostly determined by their initial 
envelope mass fraction, a behavior that dominates from 
fenv,i = 0.01 to f en v,i = 0.10 and results in a spread of 
~2 R® in radius at Gyr ages. 

We also find that the initial radius of a planet has only 
a secondary effect on its subsequent evolution; there is 
little difference after 100 Myr between a close-in mini- 
Neptune planet formed with a relatively hot start and an 
initial radius of up to 20 R® and one formed with a rel¬ 
atively cold start with an initial radius of 5 R®. Planets 
with initial radii of 10 and 20 R® converge to nearly the 
same radius within just a few Myr, thus removing one im¬ 
portant source of uncertainty in evolutionary modeling. 
The exception to this is the case of lower-mass planets 
for which mass loss dominates. 

Lower-mass (~2 M®) planets in general have larger 
radii and experience more mass loss, as well as a greater 
spread in the amount of mass loss, with representa¬ 
tive models losing 50%-80% of their envelopes, compared 
with higher-mass (4-8 M®) planets that lose 5%-20% of 
their envelopes (for representative models). In closer or¬ 
bits, low-mass planets are also much more prone to losing 
their envelopes entirely. 
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The difference in semi-major axis is relatively small, a 
factor of 2-3, between a planet that retains most of its 
envelope and a planet that loses most or all of it, but 
this distance range varies with mass. For more massive 
planets of 16 there is a transition from none of the 
envelope remaining at 2 Gyr at 0.03 AU to ~85% of it 
remaining at 0.10 AU. For low-mass planets of 2 M®, 
~3% of the envelope remains at 2 Gyr at 0.20 AU, and 
a majority of it is retained only at >0.30 AU. A simi¬ 
lar effect occurs in the radius near the inner boundaries 
of these “transition regions” as the loss of most of the 
envelope mass causes the radius to shrink rapidly. 

Increasing metallicity results in a slight increase in the 
radius of a planet predominately due to the greater opac¬ 
ity and slower cooling of the envelope. At Gyr ages, this 
difference amounts to 0.4 R® between lx and 30x so¬ 
lar metallicity. The relatively small expected range of 
metallicities for mini-Neptune atmospheres (>30x solar) 
means that metallicity is likely a small factor in varia¬ 
tions in the real planet population. 

The specific prescription for mass loss in our calcu¬ 
lations is important. It is not possible to fit most of 
the low-mass, low-density close-in mini-Neptune planets 
without mass loss, and changing the mass loss efficiency 
and the time of the onset of mass loss (the migration 
time, tM, which may be long after formation due to sec¬ 
ular processes) can result in variations in radius of several 
tenths of an Earth radius at Gyr ages, with greater varia¬ 
tions at low masses. There is significantly more variation 
in the amount of mass loss between models than in ra¬ 
dius, but the envelope mass is a relatively small fraction 
of the total mass in our models, so the variation in mass 
loss does not result in a significant change in the expected 
mass distribution of planets over time. 

In fitting observed planets in the Kepler-11 system, 
we find significantly lower envelope fractions than previ¬ 
ous studies—by about a factor of 3 for Kepler-lld and 
-lie, for example. For Kepler-lld, we find a present-day 
f env = 2.6 ± 0.3%, compared with f env = 8.4^4% com¬ 
puted bv lLopez et al.l (2012). Similarly, for Kepler-lie, 
we find f env = 6.5 ± 0.5%, compared with 17.2^'^% 
in the Lopez study. In addition to our reanalysis, we 
also achieve a preliminary fit to the latest mass and ra¬ 
dius measurements of Kepler-llg of f env = 2.6 ± 0.3%. 
Further work is needed to refine these results through 
a better understanding of the equations of state, heat 
capacity, opacities, and other factors that contribute to 
super-Earth and sub-Neptune planetary evolution. 
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However, the lower-mass and lower-density planets of 
Kepler-llb, -11c, and -Ilf are challenging to model at 
their current locations. We successfully model Kepler-llf 
with f env ~ 0.0025, an order of magnitude lower than 
the Lopez study. However, for Kepler-llb, we find no 
models that retain any H 2 -He envelope, which is required 
to account for their radii, and for Kepler-11c, we find that 
either a very low mass loss efficiency or a late migration 
time is required to retain a H 2 -He envelope. 

For Kepler-llb, the only solution we find that is con¬ 
sistent with observations is that it is a “water-world” 
with an ice core and a thick steam atmosphere and no 
hydro gen or helium (|Howe fc Burrowsll2012t iHowe et al.1 
I2014T) . While further study is needed on the evolution of 
a “water-world” in such an environment, we construct a 
plausible model to explain both Kepler-llb and Kepler- 
11 c—likely our most plausible model for both planets—in 
which the two planets formed beyond the snow line and 
then migrated inwards to their current orbits due to secu¬ 
lar processes at an age of ~300 Myr (but not significantly 
earlier). Kepler-llb, being lower in mass, either formed 
without a H 2 -He envelope or lost it entirely after migra¬ 
tion, while Kepler-llc retains a small hydrogen envelope 
on the order of f env ~ 0.007. Additional work is needed 
to investigate the possibility of extending this solution to 
all six planets in the Kepler-11 system. 

Further study is required to determine the structures 
of very low-mass, low-density planets and to more ac¬ 
curately determine the effects of the properties of the 
atmospheres (e.g. anisotropic heating, thermal redistri¬ 
bution, and clouds) on super-Earth and mini-Neptune 
planetary evolution. More accurate mass measurements 
in particular are also needed for observed planets, given 
that the deduced compositions and evolutionary histories 
can be very sensitive to mass, to accurately fit structural 
and evolutionary models to these objects. 
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